﻿#include "BAB_Lagrangian.h"
#include "assistant.h"
#include <map>
#include <numeric>
#include <iostream>
#include <cfloat>

using namespace std;
/*****************************************
 *
 * distance 的外层代表有 i 个用户，内层代表每个用户与每个设施的距离
 * fnum 记录的是需要选择的设施数量
 * ftype 存储的是设施编号与设施类型的对应关系: 'F' 后者 'C'
 *
 *****************************************/
pair<bool, vector<int>*> Lagrangian(unordered_map<int, unordered_map<int, double>>& distance, int k, int fnum, unordered_map<int, char>& ftype, vector<int>& fac_id, double alpha, unordered_set<int>& selected, vector<int>& unselected, int max_iter, double app_rate, ostream& fout, wh::reference_info_list& user_info)
{
	/* 首先进行初始化 */
	int user_size = distance.size();
	int fac_size = distance[0].size();
	double pre_beta = DBL_MIN;
	int iter_num = 1;  // 记录当前迭代的次数
	double v = 0;  // 初始化拉格朗日乘子 v
	vector<double> u(user_size, 0.0); 
	for (int i = 0; i < u.size(); i++)  // 初始化拉格朗日乘子 u_i
	{
		double min_dist = DBL_MAX;
		for (auto& p : distance[i])
		{
			if (p.second < min_dist)
				min_dist = p.second;
		}
		u[i] = min_dist;
	}

	vector<int>* res = nullptr;
	double pre_dist = alpha;
	/* 迭代 */
	do
	{
		// 首先选择 |F| 个最小的设施
		map<int, double> pj;
		{  // 代码块用于找出其中 |F| 个最小的设施，可以优化
			map<double, vector<int>> pj_dual;  // 用于 pj 中的最大元素
			for (auto i : fac_id)
			{
				// 判断 |F| 的类型
				switch (ftype[i])
				{
				case 'F':
					for (int j = 0; j < user_size; j++)
						pj[i] += min(0.0, distance[j][i] - u[j]);
					break;
				case 'C':
					for (int j = 0; j < user_size; j++)
						pj[i] += min(0.0, distance[j][i] - u[j]);
					pj[i] += v;
					break;
				}
			}

			if (unselected.size())  // 有一些点不能选择
			{
				for (auto i : unselected)
					pj.erase(i);
			}

			for (auto& p : pj)
			{
				if(selected.count(p.first) == 0)
					pj_dual[p.second].push_back(p.first);
			}

			while (pj.size() > fnum)
			{
				pj.erase(*(pj_dual.rbegin()->second.rbegin()));
				pj_dual.rbegin()->second.pop_back();
				if (pj_dual.rbegin()->second.size() == 0)
					pj_dual.erase(prev(pj_dual.end()));
			}
		}
#ifdef DEBUG
		std::cout << "\t被选出来的设施为: ";
		for (auto& p : pj)
			cout << p.first << ' ';
		cout << endl;
#endif // DEBUG


		// 确定 x_ij(用户是否向某个设施请求服务) 和 y_j(某个设施是否建立) 的值
		vector<unordered_map<int, char>> user;
		for (int i = 0; i < distance.size(); i++)  // 遍历每个用户
			user.push_back(unordered_map<int, char> ());
		for (auto& p : pj)
		{
			int f_id = p.first;
			for (int i = 0; i < distance.size(); i++)
			{
				if (distance[i][f_id] - u[i] <= 0)
					user[i][f_id] = 1;
			}
		}

		// 用于更新 u 的 theta
		vector<int> theta(user_size, 1);
		for (int i = 0; i < user_size; i++)
		{
			for (auto& p : user[i])  // user 中记录着每个用户向每个设施请求服务的情况，如果请求服务则为 1，否则不记录
				theta[i] -= p.second;
		}

		// 用于更新 v 的 eta
		int eta = k;
		for (auto& p : pj)
		{
			if (ftype[p.first] == 'C')
				eta -= 1;
		}

		// 计算 beta（对偶问题的当前值）和 alpha （原问题的解）
		double beta = 0;
		for (int i = 0; i < user_size; i++)
		{
			for (int j: fac_id) 
			{
				if (user[i].count(j))
				{
					beta += distance[i][j] - u[i];
				}
			}
		}
		fout << "距离为: " << beta;
		for (auto& p : pj)
		{
			if (ftype[p.first] == 'C')
				beta += v;
		}
		beta += accumulate(u.begin(), u.end(), 0.0);
		beta -= v * (double)k;


		if (beta > pre_beta)
			pre_beta = beta;
		if (alpha < (1+app_rate)*beta)  // 如果不满足优化的最小需求，直接返回
			return make_pair(false, res);

		/* 判断此时的解是否取得了一个更紧的上界 */
		if (eta >= 0)
		{
			double temp_dist = calculate_original_dist(distance, pj);
			if (temp_dist < pre_dist)  // 如果得到了更好的解，记录以便返回
			{
				fout << "更新原问题的解";
				pre_dist = temp_dist;
				if (res == nullptr)
					res = new vector<int>;
				else
					res->clear();

				for (auto& p : pj)
					res->push_back(p.first);
			}
			else
			{
				fout << ";新的解大于原问题的解，无法更新原问题的解";
			}
		}
		else {
			fout << ";eta小于0，不符合条件，无法更新原问题的解";
		}

		// 计算步长
		int form = 0;
		for (auto i : theta)
			form += abs(i);

		// 判断是否取得了最优解
		if (form == 0 && eta >= 0)
		{
			fout << "第" << iter_num << "次迭代，找到最优解\n";
			fout << "被选出来的设施为: ";
			for (auto& p : pj)
				fout << p.first << ' ';
			fout << endl;
			break;
		}
		double step = (alpha - beta) / form / iter_num;
		if (iter_num < 5)
			step /= 50.0;
		else if (iter_num < 20)
			step /= 25.0;
		else if (iter_num < 30)
			step /= 15.0;
		else if (iter_num < 40)
			step /= 5.0;

		if (step >= 0 && step < 0.001)
			step *= 5;

		fout << ";对偶问题的值为: " << beta << ";步长为: " << step;

		// 更新 u 和 v
		for (int i = 0; i < distance.size(); i++)
		{
			u[i] += theta[i] * step;
		}
		v = max(0.0, v - eta * step);
		fout << ";新的 v 为: " << v << endl;
		++iter_num;
	}while (iter_num < max_iter);
	return make_pair(true, res);
}

